Ecuaciones Diferenciales Ordinarias

Práctica de Cálculo, E.U.A.T., 2008

Una ecuación diferencial ordinaria (EDO) es una igualdad que involucra a una función y/o a alguna de sus derivadas. Como por ejemplo:

y'[x] == Sin[x]    

    u''[x] == 
      u[x]

Una solución de una ecuación diferencial es una función que verifica dicha igualdad. Asi las funciones - Cos[x] y Exp[x] son soluciones de las anteriores ecuaciones respectivamente. Veámoslo

-Cos'[x] == Sin[x] 
Exp''[x] == Exp[x]

Las soluciones de las ecuaciones diferenciales no son por lo general únicas, Por ejemplo, cualquier función de la forma -Cos[x] +K, con K una constante, es solución de la primera ecuación, mientras que Exp[-x] es tambien solución de la segunda ecuación:

D[-Cos[x] + K, x] == Sin[x]
D[Exp[-x], {x, 2}] == Exp[-x]

Para poder elegir entre las posibles soluciones que tiene una ecuación diferencial es necesario añadir ciertas condiciones que tienen que verificar en algunos puntos. Dependiendo  de que estas condiciones se impongan de una forma o de otra el problema completo (ecuación diferencial + condiciones puntuales) recibe distintos nombres. Un Problema de Valores Iniciales (PVI) es buscar la solución de una ecuación diferencial que en un punto tome ciertos valores impuestos. Por ejemplo si nos dicen que busquemos la solución de la ecuacion diferencial 1 que verifica que en el x=0 vale 1, entonces la única función que verifica ambas condiciones es -Cos[x]+2

D[-Cos[x] + 2, x] == Sin[x]
-Cos[0] + 2 == 1

Sin embargo si se imponen condiciones sobre la función o sus derivadas en dos puntos distintos diremos que tenemos un Problema de Valores en la Frontera.Así, si buscamos una función que verifique la ecuación diferencial 2 y que además tome los valores 0 en x=0 y en x=1, la única función que trivialmente verifica las tres cosas a la vez es la función constantemente igual a 0. En general en número de conciciones que hemos de imponer para que podamos seleccionar una solución de una ecuación diferencial coincide con el orden de la ecuación, es decir, el orden máximo de las derivadas de la función incógnita que estan involucradas en la ecuación. Así, la ecuación diferencial 1 es de orden 1 mientras que la 2 es de segundo orden.

En esta práctica veremos como resolver con Mathematica  algunas de estas ecuaciones que surgen motivadas por problemas cotidianos tanto de manera exácta como numérica (calcula una aproximación). En ciertos casos se puede obtener la  solución de la ecuación diferencial de forma exacta, es decir, su expresión analítica. Para ello Mathematica emplea la orden DSolve. Mientras que en la mayoría de los casos la solución de cada problema particular únicamente se puede calcular de forma aproximada numéricamente.  El comando en cuestión en este caso  es: NDSolve,  siendo el resultado proporcionado es una función interpolante para la solución numérica. 
Notese que ambos comandos han experimentado mejoras en las últimas versiones de Mathematica, por lo que es posible que el mismo comando a veces de solución o no dependiendo de la versión del programa empleada para calcularla.

Comandos DSolve[] y  NDSolve[]

El comando para intentar calcular de forma exacta una solución de una EDO es DSolve

DSolve[lista de EDO y Condiciones Iniciales ,  lista de funciones incognita , Variable independiente ]

Si Mathematica es capaz de calcular dicha solución nos la mostrará en pantalla, aunque hay que decir que esto no siempre es posible, por lo que hay que obtar por intentar aproximar dicha solución en unos cuantos puntos (aproximación numérica) con lo cual al menos nos podemos hacer una idea de como será la solución buscada. Para ver mas detalles

?? DSolve

NDSolve  es quizás el más importante de los comandos numéricos de Mathematica y puede resolver de forma numérica una gran cantidad de EDOs. Su sintaxis básica es

NDSolve[lista de EDO y Condiciones Iniciales ,  lista de funciones incognita , Variable independiente e intervalo]

e intenta resolver numéricamente las EDO(s) con condiciones iniciales dadas, que es preciso especificar siempre, en lista de EDOs y Condiciones Iniciales para las funciones en Funciones incognita. La variable independiente es VariableIndependiente.  Cuando este comando tiene éxito, su salida es una lista de funciones que Mathematica representa como InterpolatingFunction[ ].

El resultado de algunos comandos de Mathematica como NDSolve[ ], se muestra como "InterpolatingFunction[intervalo,listadeDatos]" que representa una función que interpola los datos dados en el intervalo especificado.

?? NDSolve

Ejemplos

Caída libre

En cinemática, la caída libre es un movimento de un cuerpo dentro del campo graviatorio terrestre. Si en este movimiento se desprecia el rozamiento del cuerpo con el aire, es decir, se estudia en el vacío. El movimiento de la caída libre es un movimiento uniformemente acelerado. Para caídas desde alturas de sólo unos pocos kilómetros o metros, la aceleración instantánea debida sólo a la gravedad es casi independiente de la masa del cuerpo, es decir, si dejamos caer un coche y una pulga, ambos cuerpos tendrán la misma aceleración, que coincide con la aceleración de la gravedad (g). Sabemos por la segunda ley de Newton que la suma de fuerzas F es igual al producto entre la masa del cuerpo y la aceleración. La ecuación de movimiento es por tanto: F=m  v' =m x''

DSolve[F == m*x''[t],  x, t]

Supongamos un cuerpo de masa=1, que inicialmente esté en reposo a 100 metros de altura y caiga en caida libre bajo el efecto de la gravdad F= -9.8,
La pregunta que nos hacemos es entonces a que altura estará cuando transcurra un segundo.  Y  a los 4 segundos???

sol1 = 
      DSolve[{-9.8 == x''[t],  
          x[0] == 100, x'[0] == 0}, x, t]

x[t_] := 100 - 4.9 (t^2)

x[1]

Plot[x[t], {t, 0, 1}];

La viga

En ingeniería y arquitectura se denomina viga a un elemento constructivo lineal que trabaja principalmente a flexión. En las vigas la longitud predomina sobre las otras dos dimensiones y suele ser horizontal.

El esfuerzo de flexión provoca tensiones de tracción y compresión, produciéndose las máximas en el cordón inferior y en el cordón superior respectivamente, las cuales se calculan relacionando el momento flector y el segundo momento de inercia. En las zonas cercanas a los apoyos se producen esfuerzos cortantes o punzonamiento. También pueden producirse tensiones por torsión, sobre todo en las vigas que forman el perímetro exterior de un forjado. Estructuralmente el comportamiento de una viga se estudia mediante un modelo de prisma.

La teoría de vigas es una parte de la resistencia de materiales que permite el cálculo de esfuerzos y deformaciones en vigas. Si bien las vigas reales son sólidos deformables, en teoría de vigas se hacen ciertas simplificaciones gracias a las que se pueden calcular aproximadamente las tensiones, desplazamientos y esfuerzos en las vigas como si fueran elementos unidimensionales.
I1=movimiento de inercia.
E1=módulo de elasticidad de Young.
P(x)=carga ejercida.
y(x)=curvatura elástica.
La ecuación que establece la relación entre la carga ejercida sobre la viga y la curva elástica es (((d^4) y)/(dx^4))=(1/(E1*I1))P(x)

E1 = 0.02;
I1 = 1;
L = 1;
P[x_] := x - L;

Plot[P[x], {x, 0, L}]

sol2 = 
      NDSolve[{y''''[x] == 1/(E1*I1)*P[x], 
          y[0] == 0, y'[0] == 0, 
          y[L] == 0, y'[L] == 0}, y, {x, 0, L}]

Plot[y[x] /. sol2, {x,  .0, L}, 
        AspectRatio -> Automatic];

xd

La braquistócrona

Consideramos dos puntos, M y N, supongamos un hilo rígido extendido entre ellos, a lo largo del cual se desliza una bola. Según sea la forma del hilo, la bola tardará más o menos tiempo en llegar del extremo M superios al extremo inferior N. Se denomina problema de la braquistócrona a encontrar la  la curva que da al tiempo mínimo de todos los posibles. La ecuación a la que se reduce a estudiar es y[x](1+y'[x]^2)=c. Se sabe que la solución de este problema es una cicloide que en se puede expresar de forma paramétrica como {c(t-Sin[t])/2,-c (1-Cos[t])/2}. Veamos entonces como se puede obtener de forma numérica.

sol3 = 
      NDSolve[{y'[x] == Sqrt[1/y[x]  - 1], 
          y[(1 - Sin[1. ])/2] == (1 - Cos[1. ])/2}, 
        y[x], {x, 0.01, Pi/2}]

Plot[-y[x] /. sol3, {x, 0.01, 
          Pi/2}];

Ejercicios

1) El perfil de un cable colgado entre dos puntos de la misma altura viene dado por la solución del Problema de Valores en la Frontera u''[x] = k (Sqrt[ 1 + u'[x]^2])  u[0]= 1, u[1]=1 con k una constante positiva que depende del peso del cable y de la tensión aplicada a dicho cable. Si suponemos que estamos empleando siempre el mismo material indica simulando con distintos valores de k que soluciones de las obtenidas se verían sometidas a una tensión mayor que otras.